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WIMP Dark Matter and the QCD Equation of State 
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Weakly Interacting Massive Particles (WIMPs) of mass m freeze out at a temperature Tj ~ m/25, 
i.e. in the range 400 MeV - 40 GeV for a particle in the typical mass range 10 - 1000 GeV. The 
WIMP relic density, which depends on the effective number of relativistic degrees of freedom at T/, 
may be measured to better than 1% by Planck, warranting comparable theoretical precision. Recent 
theoretical and experimental advances in the understanding of high temperature QCD show that 
the quark gluon plasma departs significantly from ideal behaviour up to temperatures of several 
GeV, necessitating an improvement of the cosmological equation of state over those currently used. 
We discuss how this increases the relic density by approximately 1.5 - 3.5% in benchmark mSUGRA 
models, with an uncertainly in the QCD corrections of 0.5 - 1 %. We point out what further work 
is required to achieve a theoretical accuracy comparable with the expected observational precision, 
and speculate that the effective number of degrees of freedom at T/ may become measurable in the 
foreseeable future. 
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A leading candidate for the primary constituent of 
dark matter is the Weakly Interacting Massive Particle 
(WIMP, see e.g. WIMPs of mass m generically freeze 
out at a temperature of around Tf ~ r7i/25 GeV, and so 
for typical masses in the range 10 ~ 1000 GeV, Tf lies 
between 400 MeV and 40 GeV. The equation governing 
the relic density depends on both energy and entropy 
densities, and so the WIMP relic density is sensitive to 
the equation of state of the Universe in this temperature 
range. It is normally supposed that, above the QCD 
confinement critical temperature of Tc ^ 200 MeV, the 
plasma is weakly interacting because of asymptotic free- 
dom, and can be treated as an ideal gas, and this as- 
sumption is built into at least two of the best-known 
SUSY WIMP relic density packages, DarkSUSY and 
MicrOMEGAs 0. 

However, intensive non-perturbative studies and ex- 
periments revealed that the high-temperature QCD 
plasma still departs significantly from an ideal gas at 
temperatures several orders of magnitude higher than 
Tc- It is therefore worth re-examining the equation of 
state for the Universe in this critical regime, and qual- 
itative investigations have been made before 0. Usin, 
recent progress in both the deconfined and confined 
phases, we construct an improved equation of state and 
investigate the consequences for the WIMP relic density 
in one of the dark matter packages, DarkSUSY. Dark- 
SUSY currently uses an equation of state based on 
work in Refs. 00. Replacing it with the improved ver- 
sion gives upward corrections of over 2% to the relic den- 
sities in a set of benchmark mSUGRA models ^3 j ^.nd 
larger effects for models with lower Tf. 

Although small, this correction is becoming significant 
in the new era of precision cosmology: the cold dark mat- 



ter (CDM) density is determined to better than 10% in 
single-field inflation models from a combination of the 
Cosmic Microwave Background (CMB) angular power 
spectrum and observations of the galaxy power spectrum 
by 2dF [13 and by SDSS 1^. Planck promises to do 
much better, with one estimate giving a determination 
to better than 1% [l^ . 

The calculations we use still involve simplifications and 
systematic errors, forcing us to model p{T) to some ex- 
tent, as described below. Estimating our systematic er- 
rors on the QCD equation of state to be on the order 
of 10% near Tc, we show how they propagate into un- 
certainties in the range of 0.5 - 1% in the relic density 
of WIMPs. Thus, the quest for new physics in observa- 
tional cosmology depends on further improvements in the 
accuracy of the quantitative understanding of the QCD 
equation of state. 

There has been much effort in calculating the pres- 
sure p{T) of an SU(A^c) gauge theory with Nj funda- 
mental fermions at temperature T , as reviewed in Refs. 
[l^ ll^l 0, 0| . Perturbative expansions in the coupling 
constant g of quantities such as the pressure converge 
badly, and particularly for a strongly-coupled theory like 
QCD. Strictly perturbative expansions, even when ex- 
panded to 0{g^ Ing) 0, seem to converge well only at re- 
markably high temperatures, above lO^Tc, in sharp con- 
tradiction to the ideal gas assumption usually made in 
cosmology. 

In the high temperature regime, progress has been 
made using perturbative finite-temperature dimensional 
reduction (DR) (see and references therein). By con- 
structing a sequence of effective theories for the scales 
2nT, gT and g^T, the last of which has to be treated 
non-perturbatively, one can get results for arbitrary Nf 
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of massless quarks applicable down to a few times Tc- 
By fitting for the non-perturbative and as yet unknown 
0{g^) coefficient, the authors of Ref. 0] were able to 
match their calculated pressure reasonably well to pure- 
glue lattice data near the critical temperature. 

Around the transition, there now exist lattice calcula- 
tions for the pressure and energy density for Nf = 0, 2, 3 
degenerate quark flavours, as well as first data for Nf = 
2 H- 1, i.e. two light and one heavier flavour. The pseudo- 
critical temperatures Tc{Nf) (defined as the peak of a 
susceptibility) are currently given as Tc(0) = 271 ± 2 
MeV, Tc(2) = 173 ± 8 MeV and Tc(3) = 154 ± 8 MeV 
[l^ . It has to be stressed, however, that only the pure 
glue case has been extrapolated to the continuum. Based 
on experience with that theory, the dynamical fermion re- 
sults are estimated to display a systematic error of about 
15% at the currently available lattice spacing. Correc- 
tions due to quark masses deviating from the physical 
ones appear to be negligible in comparison [l^ . 

Finally, below the phase transition the hadronic res- 
onance gas model, which treats the plasma as an ideal 
gas of mesons, baryons and their excited states, matches 
reasonably well to lattice data . 

We now review relic density calculations. Consider 
a particle of mass m and number density n, undergo- 
ing annihilations XX — > . . . with total cross-section cr, 
assumed to be a typical weak interaction cross-section, 



proportional to Gj^. Then [foj 
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where rioq is the equilibrium number density, H — a/a 
is the Hubble parameter, and WM0I is the M0ller veloc- 
ity which is a relativistic generalisation of the relative 
velocity of the annihilating particles ^3 ■ 

In order to solve the equation it is convenient to con- 
vert the time variable to cc = m/T, and to measure the 
relic abundance in terms of y — n/ s, where s is the en- 
tropy density. If the total entropy S = sa^ is conserved, 
then we can write 



cq/ 



(2) 



dY _ I ds 
dx SH dx 

This adiabaticity assumption is violated if the QCD tran- 
sition is first order, but it is most likely to be a cross-over 
transition at the low chemical pote ntials which are rele- 
vant for the early Universe mill. 

Using the Friedmann equation = SnGp/S, and 
defining effective numbers of degrees of freedom for the 
energy and entropy densities, s = {p + p)/T, through 



one finds an approximate solution 
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where Tf is the freeze-out temperature, defined to be the 
temperature at which the relic abundance is a certain 
factor (taken to be 2.5) above the equilibrium abundance. 
We see that the relic density depends on the parameter 
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It is through this parameter that the QCD equation of 
state influences the relic density. 

In an ideal gas at temperature T, a particle of mass 
rrii — XiT contributes to gcff j ^cff the amounts 
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where po ^.nd Sq are the energy and entropy densities for 
a free massless boson. 

Interactions correct the ideal gas result, and goff, hos 
have to be extracted from calculations of the energy and 
entropy, Eqs. lO. Note that in the relevant temperature 
range 40 to 0.4 GeV, 16 out the 18 bosonic degrees of free- 
dom are coloured, with the coloured fermionic degrees 
of freedom dropping from 60/78 to 36/50. The domi- 
nant corrections are therefore expected to come from the 
coloured sector of the Standard Model, weak corrections 
are moreover suppressed by the boson masses and negli- 
gible. 

The rehc density codes DarkSUSY and MicrOMEGAs 
use identical equations of state, developed in Refs. 
0, H, ■ Below Tc the hadronic degrees of freedom are 
modelled by an interacting gas of hadrons and their res- 
onances, while above Tc the quarks and gluons are taken 
to interact with a linear potential Vqir) = Kr, with a 
phenomenologically motivated value K = 0.18GeV^, de- 
rived from the slope of Regge trajectories. In this model, 
the pressure is already very close to ideal at temperatures 
above 1.6 GeV. All other Standard Model particles are 
free. 

In this work we also take the ideal gas contributions 
for the particles of the Standard Model, with masses 
given by the Particle Data Group central values [20j . 
In the confined phase quarks and gluons are replaced 
by hadronic models described below. In the deconfined 
phase, the contribution to the pressure of the coloured 
degrees of freedom is scaled by a function /(T), defined 
to be the ratio between the true QCD pressure p{T) and 
the Stefan-Boltzmann result psB for the same theory, 
/(T) = p{T)/pqb{T). This correction factor is derived 
from lattice jlj| and perturbative Q calculations for 
Nf — 0, and uses an approximate universality in the pres- 
sure curves for different Nf observed by Karsch et al. [l8j |. 
Near the transition, the lattice-derived curves for /(T) 
have the approximate form /(T, Nf) = fQCT){T /Tc{Nf)), 
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where Tc{Nf) is the critical temperature for the theory 
with Nf hght fermion flavours. Besides this, there ap- 
pears to be only negligible A^/-dependence within the 
current numerical accuracy. We are therefore motivated 
to neglect quark mass effects, take the Nf — lattice 
data (in the continuum limit) and the Nf — DR for- 
mula of Kajantie et al. [3, and scale the temperature 
dependence by Tc(3)/Tc(0). The correction factors are 
matched at 1.2 GeV using the undetermined 0{g^) pa- 
rameter, which we take to be 0.6755, close to the value 
0.7 used in Ref. At higher temperatures one crosses 
the c and b mass thresholds. However, for Nf > 3 the 
0{g^) fitting parameter is unknown, and the appropriate 
critical temperatures are also unknown. Hence we be- 
lieve that scaling the Nf = result is the best we can do 
at present. We will discuss later how improvements can 
be made. 

In the confined phase we label our model equations of 
state (EOS) A,B and C. EOS A ignores hadrons com- 
pletely, as the lattice shows that /(T) very rapidly ap- 
proaches zero below Tc. EOS B and C model hadrons 
as a gas of free mesons and baryons. It was noted in 
Ref. |6| that such a gas, including all resonances, gives 
a pressure which fits remarkably well to the Nf = 2-1-1 
lattice results although, as the authors themselves point 
out, this result should be treated with caution as the 
simulations are not at the continuum limit. We include 
all resonances listed in the Particle Data Group's table 
mass_width_02 .mc |2^. 

We make a sharp switch to the hadronic gas at a tem- 
perature Thg- For our EOS B we take Trg = = 154 
MeV, and for EOS C we take Thg = 200 MeV and 
Tc = 185.5 MeV, values chosen to give as smooth a curve 
for hcs as possible. The effects of these equations of state 
on the relic densities turn out to differ by less than 0.3% 
in the relevant temperature interval, so in the following 
we concentrate on EOS B. 

Before presenting our results we note that fic^i^ is di- 
rectly proportional to the entropy density today sq = 
{2iT^/A5)hes{T.y)T^, and that this must be determined 
as accurately as possible. The photon temperature Tj — 
2.725 ± 0.001 is very accurately measured but the 
contribution from neutrinos requires a separate freeze-out 
calculation. Recent work gives /leff(Ty) = 3.9172. 
We find that taking freeze-out temperatures to be 3.5 
MeV for and i^t and 2 MeV for Ve , as recommended in 
Ref. 0, gives /leS (Ty) = 3.9138, which is accurate enough 
for our purposes. DarkSUSY uses /leff(Ty) = 3.9139. 

In Figs. 1 we plot for our EOS B the effective num- 

1 /2 

bers of degrees of freedom hcs{T) and (T) defined 
in Eqs. Q and lO, compared with those used in Dark- 
SUSY |2] and MicrOMEGAs 0. The spike in our gl^'^{T) 
is an artefact of the matching of the scaled lattice data 
and the hadronic equation of state (i.e. the first deriva- 
tive jumps). Since for physical QCD the transition is a 
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FIG. 1: Degrees of freedom factors gV^{T) and hcs{T), de- 
fined in Eqs. |3] and |S] for our equation of state B (solid) 
compared with those currently used in DarkSUSY and Mi- 
crOMEGAs (dashed). The spike in the upper panel is a nu- 
merical artefact without physical impact, see text. 



smooth crossover, this spike is unphysical. However, it 
has no noticeable influence on the freeze-out of WIMPs 
at higher temperatures. 

In Table d we exhibit the effect of the new equation 
of state on the density of relic neutralinos Xj for the 
mSUGRA models used to test DarkSUSY in the standard 
distribution jS^]. We find changes of about 1.5-3.5%. In 
order to quantify the effect of uncertainty in the lattice 
data, we introduce two new models B2 and B3, which 
are constructed by scaling the lattice curve by 0.9 and 
1.1 respectively, and then adjusting the 0{g^) parame- 
ter in the DR pressure curve so that it meets the scaled 
lattice curve at T = 4.43Tc. Thus a 10% uncertainty 
in the lattice pressure curve translates to an uncertainty 
in the relic density in the range 0.5 - 1%. Note that 
the lowest freeze-out temperature in the table is about 
4 GeV. This corresponds to more than 20rc, where the 
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TABLE I: Relic densities in benchmark models f]3 . calcu- 
lated with DarkSUSY. Displayed: neutralino mass nix, freeze- 
out temperature Tf, flch? computed using the current Dark- 
SUSY equation of state, Qch^ using our QCD-corrected Equa- 
tions of State 65^3!, and percentage changes A. 



QCD corrections to 5, are around 5%. Evidently, WIMPs 
with freeze-out temperatures closer to Tc, such as sterile 
neutrinos Q, would be affected more strongly. 

To conclude: by updating the equation of state of 
the Standard Model in the light of recent developments 
in high temperature QCD, we have found differences in 
WIMP rehc density calculations of a few per cent for 
the benchmark models of Ref. a small but not in- 
significant effect. These models freeze out well above Tc, 
and so it is likely that there are other models exhibiting 
greater changes. In single field inflation models, the 68% 
confidence limit on the dark matter density is already 
OJi^ = 0.127 ± 0.017 from the WMAP First Year data 
[ill , while Planck has been estimated to be able to reach 
a level of An^h^ = 0.0011 In order to benefit from 
this accuracy, it has to be matched by that of theoretical 
calculations of the QCD equation of state. To obtain an 
accuracy of 1% we estimate that we require lattice data 
accurate to 10% at about T ~ 4Tc, for realistic quark 
masses, and in the continuum limit. Further progress 
is possible: the 0{g^) term is in principle calculable but 
requires a 4-loop computation in lattice perturbation the- 
ory 1^ . Moreover, the b and c quark mass thresholds may 
affect the QCD corrections. We also note that the equa- 
tion of state at a few Tc is also experimentally accessible 
in current and future heavy ion experiments, which gives 
them a direct cosmological motivation. 

Finally, given a sufficiently accurate g*, we speculate 
that with a combination of measurements of the neu- 
tralino mass and cross-sections in direct detection exper- 
iments and at LHC and the ILC, the effective number of 
degrees of freedom in the early Universe at a few GeV will 
become experimentally accessible, motivating attempts 
to improve further on this work. 
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